Analyses of swisscom data

Grid preparation

Ancillary data

Municipality boundaries

st_gg21 <- st_read("data-raw/swisstopo/swissboundaries3d_2021-07_2056_5728/SHAPEFILE_LV95_LN02/swissBOUNDARIES3D_1_3_TLM_HOHEITSGEBIET.shp",
                   as_tibble = TRUE) %>% 
  st_zm(drop = TRUE, what = "ZM") %>% 
  # move to LV03
  st_transform(21781) %>%
  # filter(! BFS_NUMMER %in% c(2391, 5391, 5394)) %>% 
  # exclude lakes
  # filter(OBJEKTART != "Kantonsgebiet") %>% 
  # exclude FL & enclaves
  filter(ICC == "CH") %>% 
  select(BFS_NUMMER, NAME, KANTONSNUM, GEM_TEIL) %>% 
  rename(GMDNR = BFS_NUMMER,
         GMDNAME = NAME,
         KTNR = KANTONSNUM) %>% 
  arrange(GMDNR)

write_rds(st_gg21, "data/swisstopo/st_gg21.Rds")
st_gg21 <- read_rds("data/swisstopo/st_gg21.Rds")

STATPOP

Could be used to define offset for st_make_grid

statpop_20 <- read_delim("data-raw/BfS/ag-b-00.03-vz2020statpop/STATPOP2020.zip", 
                         delim = ";", escape_double = FALSE, trim_ws = TRUE)[] %>% 
  mutate_all(as.integer)

# statpop_20 %>% slice(1:10) %>% View()

write_rds(statpop_20, "data/BfS/statpop_20.Rds")
offset <- read_rds("data/BfS/statpop_20.Rds") %>% 
  summarise(minX = min(X_KOORD),
            minY = min(Y_KOORD))

Grid

Generate whole country

100 m grid, with lower left corner defined using STATPOP cells.

Also, adding ID for each cell (based on row number).

country <- st_gg21 %>% 
  st_make_grid(cellsize = 100, 
               offset = c(offset$minX, offset$minY),
               square = TRUE) %>% 
  st_sf() %>%
  st_cast("POLYGON") %>% 
  mutate(ID1 = row_number()) %>% 
  relocate(ID1)

# st_write(country, "data/grid/country.gpkg")
write_rds(country, "data/grid/country.Rds")
country <- read_rds("data/grid/country.Rds")

Generate cities

Clipped using swisstopo municipality data.

Added second city specific ID and area.

Zuirich

zurich <- country %>% 
  st_intersection(st_gg21 %>% filter(GMDNR == 261))  %>%
  mutate(ID2 = row_number(),
         area = st_area(.)) %>% 
  relocate(ID1, ID2, area)

Geneva

geneva <- country %>% 
  st_intersection(st_gg21 %>% filter(GMDNR == 6621))  %>%
  mutate(ID2 = row_number(),
         area = st_area(.)) %>% 
  relocate(ID1, ID2, area)

Basel

basel <- country %>% 
  st_intersection(st_gg21 %>% filter(GMDNR == 2701))  %>%
  mutate(ID2 = row_number(),
         area = st_area(.)) %>% 
  relocate(ID1, ID2, area)

Lausanne

lausanne <- country %>% 
  st_intersection(st_gg21 %>% filter(GMDNR == 5586))  %>%
  mutate(ID2 = row_number(),
         area = st_area(.)) %>% 
  relocate(ID1, ID2, area)

Bern

bern <- country %>% 
  st_intersection(st_gg21 %>% filter(GMDNR == 351))  %>%
  mutate(ID2 = row_number(),
         area = st_area(.)) %>% 
  relocate(ID1, ID2, area)
# st_write(bern, "data/grid/bern.gpkg")
write_rds(bern, "data/grid/bern.Rds")

Example

Municipality Bern coverage:

bern %>% 
  st_geometry() %>% 
  qtm(borders = "red", fill = NULL)

There are few examples of grid cells with very small area.
These could be excluded if no population was found there?

bern %>% 
  st_drop_geometry() %>% 
  arrange(area) %>% 
  slice(1:10)
       ID1  ID2            area GMDNR GMDNAME KTNR GEM_TEIL
1  4239626  401 0.1281652 [m^2]   351    Bern    2        0
2  4403056 5188 0.1679092 [m^2]   351    Bern    2        0
3  4312544 2876 0.2880109 [m^2]   351    Bern    2        0
4  4319505 3098 0.3169968 [m^2]   351    Bern    2        0
5  4253464  702 0.6640466 [m^2]   351    Bern    2        0
6  4225622  139 1.0652267 [m^2]   351    Bern    2        0
7  4225603  123 1.0984176 [m^2]   351    Bern    2        0
8  4399516 5082 1.9663071 [m^2]   351    Bern    2        0
9  4270930 1374 2.5263870 [m^2]   351    Bern    2        0
10 4253513  735 3.3474380 [m^2]   351    Bern    2        0

Area estimates

nrow(bern)
[1] 5487
nrow(zurich) +
  nrow(geneva) +
  nrow(basel) +
  nrow(lausanne) +
  nrow(bern)
[1] 24116